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Abstract. There are successful approaches to explain the formation of the tachocline by a poloidal magnetic 
field in the solar core. We present here the first MHD simulations of the solar tachocline which self-consistently 
include the meridional circulation. We show that the meridional flow significantly changes the shape and the 
characteristics of the tachocline. We find that after the inclusion of the meridional circulation, a tachocline can be 
formed even when the poloidal field lines are crossing the boundary between the radiative zone and the convection 
zone. We also discuss the effects of the magnetic Prandtl number as well as of the magnetic Reynolds number on 
the properties of the tachocline. The tachocline is much thinner for higher magnetic Reynolds numbers and/or 
lower magnetic Prandtl numbers. We expect that a poloidal magnetic seed field of around 1 G will be sufficient to 
produce the tachocline of the Sun. However, the model requires the initial magnetic field to be in a narrow range 
for satisfying tachocline solutions. The simulations including a stable temperature gradient produce a shallower 
as well as slower meridional circulation than the ones without it, as desired from the Lithium abundance at the 
solar surface. 
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1. Introduction 

The term tachocline refers to a very thin layer at the 
base of the convection zone, which separates the latitu- 
dinal differential rotation inside the convection zone from 
the uniformly rotating radiative zone. Latest observations 
put the position of its center at around 0.69i?Q and esti- 
mate its thickness to be slightly less than 0.05i?Q. There is 
also general belief that the tachocline is prolate and is also 
thicker near the pole. However, the observations for higher 
latitudes have large error bars, and hence the evidence for 
the latitudinal dependence of the tachocline is not yet con- 
clusive (Schou et al. 1998, Antia et al. 1998, Charbonneau 
et al. 1999, Basu & Antia 2001). The facts were verified 
qualitatively in frequency- versus-radius plots presented by 
e.g. Eff-Darwich et al. (2002) and Howe (2003). 

There have been many attempts to explain the for- 
mation of such a thin shear layer inside the Sun. 



Spiegel & Zahn (1992) were the first to give a purely hy- 
drodynamic model for the tachocline. |Elliott (1997)| ex- 
tended this analysis numerically also including meridional 
circulation. These models were successful in explaining 
the stability of the tachocline using anisotropic turbulence 
but the equilibrium thickness of the tachocline obtained 
was much thicker than latest observational estimates. Also 
later works such as that of Miesch (2003) confirm that a 
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purely hydrodynamic model is insufficient to explain the 
tachocline. 

Riidiger & Kitchatinov (1997) studied the tachocline 
under the influence of a magnetic field. They be- 
lieved that the existence of a weak seed field of 
about 0.01 to 1 mG in the solar core is highly 
likely. And this weak internal magnetic field will be 
mainly responsible for the formation of the tachocline. 
Several other authors like Gough & Mclntyre (1998) 



MacGregor fc Charbonneau (1999)[ and |Garaud (2002) 
have used this idea of a weak poloidal seed field forming 
the tachocline. 

Here we present a magnetohydrodynamic model, 
which self-consistently calculates the meridional flow. In 
this work, we show how the inclusion of the meridional 
circulation changes the shape of the tachocline and dis- 
cuss the effects of varying magnetic Prandtl number and 
magnetic Reynolds number on the structure and the prop- 
erties of the tachocline. Finally, we discuss the effect of the 
stable temperature gradient on the structure of the merid- 
ional flow. 



For a very simple magnetic shear flow model, the rela- 
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has been derived by Riidiger & Kitchatinov (1997) for 
the fractional thickness of any magnetically induced 
tachocline layer with the Hartmann number 



Ha 



BpR 



(2) 



implying that rather thin tachoclincs of say 0.03Rq al- 
ways require Ha ~ 10 3 . The Hartmann number also ap- 
pears in the relation between the induced toroidal mag- 
netic field B$ and the given poloidal field amplitude Bq, 
i.e. 



So" 



QqR f Popv \ 
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(the amplification factor), with the global magnetic 
Reynolds number Rm = ^IqR 2 jr\. 

With r\ ~ 10 3 cm 2 /s as a rough estimate for the solar 
diffusivity, one finds Rm ~ 10 12 for the solar transition 
region between the convection zone and the radiative in- 
terior so that a very strong amplification of nine orders of 
magnitude should occur. The code used below only reaches 
Rm- values of order 10 5 , thus being limited on the ampli- 
fication factors. 



Forgacs-Dajka & Petrovay (2002) and Petrovay (2003) 



suggested that the tachocline, being just below the base 
of the convection zone, is more likely to be turbulent 
and hence the values of the diffusivities in this region 
are more likely to be as high as those in the convection 
zone. Studies by several authors concentrated on the sta- 
bility of the tachocline being subject to latitudinal shear 
(Watson 1981, Oilman & Fox 1997, Garaud 2001, Cally 
2003, Dikpati et al. 2003). These ideas, however, are be- 
yond the scope of the present discussion. 

2. The model 

The principle aim of this work will be to reproduce the 
observed rotation profile in the tachocline region, along 
with the nearly uniform rotation in the solar radiative 
zone below. 

In order to be consistent with the observations, we 
require that (i) the thickness of the tachocline should be 
less than 0.05-Rq, (ii) the core should settle near uniform 
rotation, (iii) the meridional circulation should be slow 
and/or should not penetrate the radiative zone deeply, 
and (iv) the internal seed field should be weak. 

We assume that the tachocline lies completely inside 
the radiative zone. The equations take axisymmetric, in- 
compressible Boussinesq form. Though we would like to 
include density stratification for more realistic simula- 
tions, our test runs with even moderate density gradients 
proved computationally too expensive. We take the outer 
radius of the computational domain as i? ut = 0.75Rq 
and the inner radius as R m — O.IRq. Due to numerical 
constraints, we restrict ourselves to spherical shells instead 
of a complete sphere. The size of the remaining inner hole 
was found to have negligible effect on the results. 



We also assume that the rotation profile in the convec- 
tion zone is independent of the dynamics in the radiative 
zone and can be prescribed as a boundary condition. We 
further assume that any physical phenomenon occurring 
in the convection zone, including a dynamo, does not bear 
an effect on the dynamics within the radiative zone and 
the tachocline. For the first set of simulations, the effect 
of the temperature gradient and hence buoyancy force on 
the fluid is also neglected. We will present the equations 
and results involving the temperature in Section 5. 

The magnetohydrodynamic equations employed are 

du 1 

— = -(«• V)u- VP + i/AmH (V x B) x B, (4) 

at fioP 



dB 4 
dt 



= [V x (u x B)]f - tj[V x (V x B)]„ 



(5) 



with V • u = and V ■ B = where v and r\ are the 
constant viscosity and magnetic diffusivity, respectively. 
Other symbols have their usual meaning. The poloidal 
magnetic field is maintained time- invariant. This approxi- 
mation is valid as in the Sun, due to the small 77, the mag- 
netic diffusion time is believed to be several times larger 
than the present solar age, whereas the stable tachocline 
solution is reached in our simulations in much shorter 
times. 

The equations are solved in a spherical shell geometry 
usin g the spectral cod e developed by Hollerbach (1994) 
and Hollerbach (2000) The equations are scaled in time 



by the magnetic diffusion time, Tdis — i? 2 ut /?7, in velocity 
by rj/R out , and the magnetic field is normalized by Bo, 
thus generating the dimensionless Lundquist number 



So 



B R o 



(6) 



This reduces the equations to very few free parameters, 
namely r out = R out /RQ, r in = R in /R Q , Rm, Pm = v/r) 
and So- 

Given the very high turbulent viscosity in the convec- 
tion zone, we can safely assume that the convection zone 
rotation profile acts like a rigid outer boundary on the 
tachocline. Therefore, the rotation profile of the angular 
velocity at the outer boundary is chosen to match the ob- 
served rotation profile in the bulk of the convection zone, 
which is also used by MacGregor & Charbonneau (1999) 
and is maintained as a rigid boundary condition. It is given 

by 



fl out = n (l - 0.1264cos 2 (9 - 0.1591cos 4 



(7) 



where 9 is the co-latitude 1 . The remaining boundary con- 
ditions on the velocity are maintained stress-free. For the 
magnetic field, we set vacuum boundary conditions at the 
inner as well as the outer boundary. 



1 |MacGregor fc Charbonneau ( 1999) | give a positive sign for 
the cos term which appears to be a typo. The negative sign 
used here is essential for a solar-like differential rotation. 
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Due to computational limitations, we will use values 
of Rm not exceeding 10 5 . We assume that the magnetic 
diffusivity 77 and the viscosity v take their molecular values 
in the tachocline region, but our numerically constrained 
choice of Rm implies very high values of rj and v in the 
simulations. 

The computations evolve only the axisymmetric modes 
of the originally three-dimensional spectral code. The runs 
typically employ 60 radial Chebyshev (k) modes and 60 
latitudinal Legendre (I) modes. It has been verified by 
additional computations that the results change very little 
at different resolutions. The Chebyshev polynomials are 
ideal to resolve the radial boundary layers very well. 

The magnetic field structure used for the internal seed 
field is dipolar with no initial toroidal field, i.e. 

\r 2 sin8 88' rsin(9 Or ' J K > 

with the generating function 

A = S r 2 ( 1 - —) (l-—) sin 2 9 (9) 
\ r ou t/ v r > 

The function involves the second bracket in order to avoid 
the magnetic field lines crossing the inner boundary. This 
will introduce small curls of the poloidal magnetic field 
near the inner boundary. We eliminate the curl of the 
poloidal magnetic field, which is dominant only at the in- 
ner boundary, by equating it to zero at every time step. A 
comparison with computations including the full Lorentz 
force did not show significantly different results though. 
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Fig. 1. Results for the simulations excluding (left) and 
including (right) the meridional circulation with the mag- 
netic field being confined to the core. The solid lines are 
the iso-rotation curves, whereas the dashed and the dotted 
lines represent the contours of the positive and the neg- 
ative toroidal field strength respectively. The dot-dashed 
lines indicate the boundaries of the simulation domain. 
The time-invariant poloidal magnetic field is not shown. 
S = 1100, Pm = 1 and Rm = 10 4 . 



3. The effect of the meridional flow 

Earlier investigations have shown that even a weak 
poloidal seed field is able to produce a solar-like tachocline. 
It was found that the region near the rotation axis (i.e. 
poles) is least affected by the field and the tachocline is 
thickest in that part. Further, when the poloidal field am- 
plitude is very small, the magnetic field is unable to alter 
the rotation, and the resulting rotation profile looks very 
similar to the non-magnetic case. At the other extreme, 
if the magnetic field strength is very high, the contour 
lines of constant follow the poloidal magnetic field lines 
throughout the interior, in accordance with the theorem 
of Ferraro (1937). A solar-like tachocline is thus impos- 
sible (Garaud 2002) in either case. Our first simulations 
neglecting the meridional flow concur with most of these 
key results. 

In the simulations of this Section, we use Rm = 10 4 
and Pm = 1. In all our simulations, we found that the 
system converges to an equilibrium solution in less than 
O.lTdig. Unless mentioned otherwise, the figures and nu- 
merical values throughout this Paper refer to the final 
equilibrium solutions. The results are best represented by 
the comparison of the following two graphs. While the 
left panel of Fig. is just a reproduction of the earlier 
MacGregor & Charbonneau (1999) result, where the mag- 
netic fields are confined to the simulation domain, the 



right panel of Fig. ^ shows the same model including the 
meridional circulation. The solid lines are the iso-rotation 
curves whereas the dashed and the dotted lines represent 
the contours of the positive and negative toroidal field 
strength, respectively. While in the left graph the region 
along the rotation axis is not affected at all, we see that the 
inclusion of the meridional circulation changes the picture 
completely. 

We observe that the entire core, including the region 
near the rotation axis, has achieved nearly uniform rota- 
tion. The tachocline is formed near the outer boundary. 
In contrast with the results ignoring the meridional flow, 
the tachocline is now thinnest near the pole. In the re- 
gion near the equator, where the magnetic field influence 
is smaller, we find that the iso-rotation curves tend to be 
similar to the characteristic Taylor-Proudman flow. We 
further observed that the toroidal magnetic field strength 
is only 30% of the poloidal magnetic field strength. 

Following MacGregor & Charbonneau (1999), we also 
studied the second case, where the poloidal magnetic field 
is no more confined to the computational domain and 
the magnetic field lines are crossing the outer boundary 
of the simulation domain. While this case produced no 
real tachocline for the MacGregor & Charbonneau (1999) 
model, this is no longer a problem after the inclusion of 
the meridional flow (see Fig. [2J) ■ The initial and final states 
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Fig. 2. Results for the simulations excluding (left) and in- 
cluding (right) the meridional circulation when the mag- 
netic field lines are crossing the outer boundary. 

in terms of fractional f2 are shown as functions over radius 
for various latitudes in Fig. |5J 

We performed some test simulations with a decaying 
poloidal magnetic field. During these simulations, we ob- 
served that the magnetic field does form a tachocline dur- 
ing the period up to ~ 0.05rdig, but as the magnetic field 
decays gradually due to the high value of i] assumed in 
the simulations, the flows readjust themselves, and the 
resulting pattern will be the same as the purely hydro- 
dynamic case, i.e. it will approach the Taylor-Proudman 
flow. Fig. 01 shows the fractional f2 in an intermediate 
stage at t = 0.05rdiff for a decaying poloidal magnetic 
field, keeping other parameters the same as for the run in 
Fig. |3 The effect of forming a tachocline structure by a 
time-dependent poloidal field is thus very similar to the 
non-decaying field, but diffusivity reduces the field unnat- 
urally early, and the final tachocline state is not achieved. 
Fig. [S] shows the evolution of the poloidal magnetic field 
at three different times. 

We conclude here that the magnetic field is not only 
important for the formation of the tachocline but is also 
important for maintaining it, at least in the case of the 
high magnetic diffusivity used in these simulations. 

4. Effect of varying magnetic Prandtl number and 
magnetic Reynolds number 

4.1. Varying the magnetic Prandtl number 

Following Kippenhahn & Weigert (1994) and Stix & 
Skaley (1990), the solar value for the magnetic Prandtl 
number is of the order of 1CP 3 . Although this value was 
not achieved, the simulations were performed for various 
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Fig. 3. Fractional fl vs. fractional radius for the initial 
configuration (left) and the equilibrium solution when 
Bp i is constant (right). Different lines represent Q at dif- 
ferent co-latitudes. 
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Fig. 4. Fractional fl vs. fractional radius for a snapshot 
at 0.05rdiff (left) and 0.5rdig (right) with decaying B po \. 
Again, different lines represent at different co-latitudes. 




Fig. 5. Snapshots of B po \ in case of a decaying poloidal 
field at (from left to right) t = 0, 0.05rdiff, and 0.5rdiff. 
For every snapshot, length of arrows is proportional to 
the amplitude of i? po i at each point. The length of arrows 
is renormalized in each plot by corresponding maximum 
field amplitude. Only every fifth arrow in and r is plotted 
in order to avoid congestion. 

values of Pm in the range of 0.05-1 whereas Rm = 10 4 in 
all the cases. 

In these simulations, we observed that the horizontal 
component of the meridional flow (u$) generated is about 
2% of Ucf, at higher latitudes and is about 1% for lower 
latitudes. Similarly, the vertical component of the merid- 
ional flow (u r ) is about 0.5% (see Fig.01eft). These ratios 
are nearly independent of the value of Pm. The radial 
component of the velocity remains much smaller than the 
latitudinal component especially for lower values of Pm, 
as can be seen in Fig.|H| The Figure also shows that strong 
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Fig. 6. Results for the simulations with Pm = 1 (left) 
and 0.05 (right) showing a small section of the simulation 
domain with the components of the meridional circulation 
plotted. The arrow lengths are assigned in the same way 
as Fig. [S] Only alternate points in the radial direction are 
plotted. 
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Fig. 7. The amplitudes of the components of the merid- 
ional flow at different co-latitudes (0), normalized to the 
rotational velocity at the equator, with varying Pm (left) 
and Rm (right). 



latitudinal flows are expelled closer to the outer boundary 
for lower values of Pm. The tachocline, thus formed, is 
thinner 2 at the pole than at the equator. The thickness 
of the tachocline at the equator reduces marginally from 
O.O3Oi? to O.O26i? when Pm is changed from 1 to 0.05, 
whereas the thickness of the tachocline at the pole de- 
creases considerably. We also observe that the amplitude 
of the toroidal magnetic field remains smaller than that of 
the poloidal field and is almost independent of the choice 
of Pm. 

4.2. Varying the magnetic Reynolds number 

Now we will use a constant Pm = 1 and vary Rm. As we 
know that the solar value of Rm is about 10 12 , we ran 
simulations for Rm larger than Rm = 10 4 . The highest 
value of Rm in our simulations was 10 5 . 

In this set of simulations we again notice that the rela- 
tive amplitudes of the meridional flow are nearly indepen- 
dent of the value of Rm as shown in Fig. (right). The 

2 The thickness is defined as the distance between the outer 
boundary and the radius at which the rotation rate deviates 
by 1% from the rotation rate deep inside the core. 
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Fig. 8. Variation of So with varying Pm as a log-log plot, 
for an approximate tachocline thickness of 0.025i?Q. The 
dotted line represents a fixed Ha = 1500. 



thickness of the tachocline at the equator reduces from 
0.030i?Q to 0.015i?Q when Rm is changed from 10 4 to 
10 5 , whereas the thickness of the tachocline at the pole 
goes down from 0.012i?Q to 0.003i?Q. Correspondingly, 
the toroidal field (B^) in the tachocline region increases 
from 0.3B to l.3B . 

4.3. Effect on the Lundquist number 

It will be worthwhile to give some emphasis on the choice 
of the poloidal field strength. In our simulations, we find 
that the choice of the amplitude of the seed field value 
is very critical. A small deviation in either direction from 
the magnetic field which produces a solar-like tachocline 
either makes the iso-rotation curves similar to the non- 
magnetic case or the curves will comply with Ferraro's 
theorem, unable to produce a tachocline in both cases. At 
lower values of Pm, we require smaller seed field as shown 
in Fig. |H1 At higher values of Rm, the required value of So 
is higher (Fig. EJ. 

An empirical law governing the variation of So as a 
combined function of Pm and Rm can be deduced from 
these plots as 



S ~ 10 Rm a5 Pm' 



0.25 



(10) 



This scaling is the main result of the computations. It 
means that 



B () 



(11) 



which, of course, has the correct dimension of a veloc- 
ity. The old estimate without meridional flow Rudiger & 
Kitchatinov (1997) led to the quite different expression 



Sule et al.: A numerical MHD model for the solar tachocline 



5000 
4000 

3000 
2000 



1000 



10 



1 0" 



Rn 



Fig. 9. Variation of So with varying Rm as a log-log plot, 
for an approximate tachocline thickness of 0.028-Rq. 



Bq/ \j\i$p ~ 10 3 y/ufj/R which is a very small value due to 
the appearance of the radius R. Eqn. l(TT|l yields 



Ba 



0.1 cm/s, 



(12) 



so that a maximum field amplitude of 1 Gauss results, for 
an average density of 10 g/cm 3 . This is a much larger value 
than the mGauss values for models without meridional 
flow, but it is not an unrealistic number. By contrast with 
the old model, the toroidal field belts now have the same 
order of magnitude as the poloidal fields. 

5. Effect of a temperature gradient 

In the simulations discussed in the previous sections, we 
noted that the amplitude of the meridional circulation was 
nearly independent of the variation of Rm as well as Pm. 
In the lower latitudinal belt, where the solar dynamo is 
likely to be located, the amplitude of horizontal velocity 
ug was always around 1% of u,p, and the flow reached very 
deep layers of the shell. The Lithium abundance, how- 
ever, as observed at the solar surface suggests that the 
meridional circulation should either be very shallow or 
very slow. Otherwise Lithium would be destroyed in its 
fusion zone below 0.68R@. 

We therefore include a given temperature gradient in 
our model. We introduce a temperature fluctuation on 
top of this temperature profile. This fluctuation induces a 
buoyancy force. Hence Eqn. (0} will be modified to read 



- = -(u.V } u- 

+gp, 



VP + vl\u 



— (V x B)xB 



(13) 



where g is the vector of gravitational acceleration. In the 
context of the Boussinesq approximation and after rescal- 




0.0 0.2 



Fig. 10. Results for the simulations excluding (left) and 
including (right) buoyancy. The solid and the dotted lines 
are clockwise and anti-clockwise meridional circulation, 
respectively. Ra = — 2 x 10 7 , Rm = 10 4 . 



ing as in Section 2, we obtain a buoyancy force as Ra0r 
with the modified Rayleigh number 



Ra = 



9<xATRl u 
rj 2 



(14) 



where AT = T; n — T out and a is the coefficient of vol- 
ume expansion. The equilibrium solution obtained from 
the temperature equation alone is chosen as the mean tem- 
perature profile (normalized with its bottom value) and is 
given by 



To 



V r 



1 



(15) 



and the non-dimensional energy equation has the form 



dt 



Pm 

~p7 



Ae-M-V(6 + T ), 



(16) 



where Pr = v/x is the Prandtl number using the thermal 
diffusivity x- 

When the model is evolved from the initial state in- 
cluding the temperature equation, it takes much longer 
(in terms of magnetic diffusion times) to achieve a steady 
state solution than the simpler cases in the previous sec- 
tions. To save the computational resources, we fed the 
steady state solutions obtained in the previous sections as 
the initial condition for the simulations involving temper- 
atures. It was verified that the solutions obtained in this 
manner are identical to the solutions obtained from same 
buoyancy runs but starting with the very initial conditions 
of Section 2. 

We used Rm = 10 4 and Pm = 1 for these sets of sim- 
ulations. As the simulation domain is of radiative nature, 
we use negative values of Ra. We performed the simula- 
tions with various values of Ra. The results are presented 
in Fig. For the regime when —Ra/Rm 2 < 0.01, the 
buoyancy is unable to produce any major change in any 
of the velocity components in magnitude or in structure. 
For more negative Ra, we see a gradual change in the 
structure of the meridional circulation 3 . The circulation 



Note that -Ra/Rm 2 ~ g/(tl 2 R) 
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Fig. 11. The same as in Fig.^but with the effect of buoy- 
ancy, Ra = -2 x 10 7 , Rm = 10 4 . 

is then much shallower for higher latitudes (low 9) , as de- 
sired to explain weak mixing into the interior. The mag- 
nitude of the flow also decreases but the change is not 
drastic; see Fig. El For the lower latitudinal belt (high 
9), which is an important region for the solar dynamo, 
the decrease is stronger in a relative way and the depth of 
the circulation is also clearly reduced. We get a marginal 
increase in the amplification of the toroidal field, probably 
because the field is not advected through the entire com- 
putational domain anymore. The structure of the toroidal 
field is thus changed as well, and it is shifted towards the 
outer parts of the shell. On the other hand, the rotation 
rate at higher latitudes, at large depth in the core, be- 
comes slightly slower than that at the equator. But even 
this change is marginal and within the observational lim- 
its (|1 - fipoie/Oequatorl < 3 % for r < 0.65Rq). Angular 
velocity and toroidal field belts are shown in Fig. ^] We 
expect that even more negative Ra will further reduce the 
depth and amplitude of the meridional circulation. 

6. Summary 

In this work, we have presented the first MHD model for 
the solar tachocline which self-consistently calculates the 
meridional circulation. The consideration of the merid- 
ional flow changes the shape, the structure and the char- 
acteristics of the tachocline radically. Hence, the merid- 
ional circulation cannot be neglected while modeling the 
solar tachocline. The thickness of the tachocline in the 
outer boundary layer near the equator will be determined 
by the gradient of the magnetic field near the boundary. 
Our simulations show that the tachocline is thinner near 
the pole. But the simulation domain has uniform values of 
viscosity and magnetic diffusivity, and it will only be fair 
to say that we simulate that part of the tachocline, which 
is inside the radiative zone. Hence, a definitive conclusion 
about the thickness of the tachocline in the polar region 
cannot be drawn from these simulations. 

We further report that the tachocline is thinner at 
lower values of Pm as well as at higher values of Rm. The 
meridional circulation is nearly independent of the vari- 
ation of Rm and Pm. The amplification of the toroidal 
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Fig. 12. Plot of the amplitudes of the components of the 
meridional flow at different co-latitudes 9, normalized to 
the rotational velocity at the equator, with varying Ra. 
The dashed lines represent the corresponding values with- 
out buoyancy effect. Rm = 10 4 , Pm = 1. 



magnetic field is naturally larger at higher values of Rm, 
whereas from it is clear that the poloidal magnetic 
field amplitude required to produce a solar-like tachocline 
goes down with decreasing 77, i.e. increasing Rm. The mag- 
netic seed field required to produce a solar-like tachocline 
is a function of Rm as well as Pm. The value of the seed 
field is expected to be around 1 Gauss in the Sun, for an 
average density of 10 g/cm 3 . Again a significant change is 
noted from the simulations without meridional flow where 
even a sub-mGauss field was enough to produce a solar- 
like tachocline. The scaling of the magnetic field as a 
function of the rotation rate, v, and r\ is given in 
The toroidal magnetic field is expected to be a few orders 
higher than the poloidal seed field in the case of the Sun. 

When a stable temperature gradient is introduced 
across the shell, it makes the meridional circulation shal- 
lower as well as weaker for stronger stabilization (more 
negative Rayleigh number). This can prevent Lithium 
from reaching its fusion zone which starts just below 
tachocline. The stabilizing effect is in line with the rel- 
atively high Lithium abundance observed at the surface 
of the Sun. The toroidal magnetic field in this case is lim- 
ited to belt in the outer parts of the radiative zone. 

The relations describing the variation of various pa- 
rameters such as the Lundquist number required to form 
a solar-like tachocline, the amplification of the toroidal 
magnetic field, the amplitude of the meridional circula- 
tion etc. are based on the simulation results for a limited 
range of Rm and Pm. We hope to improve the code and 
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verify the correctness of these relations closer to the solar 
values of Rm and Pm in the near future. 

Acknowledgements. L. Kitchatinov is acknowledged for the 
discussions about the simulation results and H.M. Antia for 
clarifying the observational results. J.-P. Zahn is acknowledged 
for going through the manuscript and giving valuable sugges- 
tions. 

References 

Antia, H. M., Basu, S., & Chitre, S. M. 1998, MNRAS, 298, 
543 

Basu, S., & Antia, H.M. 2001, MNRAS, 324, 498 
Cally, P. S. 2003, MNRAS, 339, 957 

Charbonneau, P., Christensen-Dalsgaard, J., Henning, R. et al. 

1999, ApJ, 527, 445 
Dikpati, M., Oilman, P., & Rempel, M. 2003, ApJ, 596, 680 
Eff-Darwich, A., Korzennik, S. G., & Jimenez- Reyes, S. J. 2002, 

ApJ, 573, 680 
Elliott, J. R. 1997, A&A, 327, 1222 
Ferraro, V.C.A. 1937, MNRAS, 97, 458 

Forgacs-Dajka, E., & Petrovay, K. 2002, A&A, 389, 629-640 

Garaud, P. 2001, MNRAS, 324, 68 

Garaud, P. 2002, MNRAS, 329, 1 

Gilman, P., & Fox, P. 1997, ApJ, 484, 439 

Gough, D.O., & Mclntyre, M.E. 1998, Nature, 394, 755 

Hollerbach, R. 1994, Phys. Fluids 6, 2540 

Hollerbach, R. 2000, Int. J. Numer. Meth. Fluids, 32, 773 

Howe, R. 2003, Proc. SOHO 12/ GONG+ 2002, 81 

Kippenhahn, R., & Weigert, A. 1994, Stellar Structure and 

Evolution. Springer 
MacGregor, K.B., & Charbonneau, P. 1999, ApJ, 519, 911 
Miesch, M. 2003, ApJ, 586, 663 
Petrovay, K. 2003, Sol. Phys., 215, 17 

Riidiger, G., & Kitchatinov, L. 1997, Astron. Nachr., 318, 273 

Schou, J., Antia, H. M., Basu, S. et al. 1998, ApJ, 505, 390 

Spiegel, E.A., & Zahn, J.-P. 1992, A&A, 265, 106 

Stix, M., & Skaley, D. 1990, A&A, 232, 234 

Watson, M., Geophys. Astrophys. Fluid Dyn., 161, 285 



